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Abstract 

We extend the result of Yan (2012) to the rhomboidal problem with masses 1, m, 
1, m for some < m < 1. A symmetry -reduction technique presented by Roberts 
(2007) enables us to compute the linear stability of the orbit by computing a single 
number. The orbit is shown to be linearly stable except for a small neighborhood 
of m = 0. Additionally, analysis of a suitable Poincare section suggests numeri- 
cally that the orbit is nonlinearly stable as well except for small m. 

Keywords: n-body problem, binary collision, regularization, linear stability, 
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1. Introduction 

The Newtonian n-body problem involves determining the motion of n point 
masses in space. The bodies are mutually attracted by the gravitational law demon- 
strated by Newton. Specifically, if {q\,qi, ...,q n } represent the positions of the 
bodies in R k (k = 1,2, or 3) with masses {mi,m2, ...,m„} respectively, then their 
motion is governed by the system of differential equations 

where the dot represents the derivative with respect to time. 

Periodic orbits of the n-body problem are of particular interest. The most note- 
worthy recent example is the three-body figure-eight orbit discovered numerically 
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by Moore in [1], and whose existence was proven analytically by Chenciner and 
Montgomery in [2]. This orbit features three equal masses which all travel on the 
same curve in a periodic fashion. This orbit was shown to be linearly stable by 
Roberts in [3] by using the twelve-fold symmetry present in the orbit, reducing 
the question of stability to a simple numerical calculation. Robert's technique has 
since been applied to other symmetric periodic orbits in various n-body problems 
(for instance [4], [5], [6], [7], and [8]). 

Another important aspect of the n-body problem involves the study of sin- 
gularities. A collision singularity is said to occur when qi = qj for some i ^ j. 
Mathematically, this yields a zero denominator in (1). Certain classes of colli- 
sion singularities are known to be regularizable, meaning that under a suitable 
change of space and time variables the orbit can be analytically continued past the 
collision time. In particular, binary collisions (where exactly two bodies collide) 
can be regularized. It is also possible to regularize simultaneous binary colli- 
sions (where multiple pairs of bodies collide at the same time) in certain cases 
when symmetries are present in the orbit. In [9] and [10], Sweatman analyzes a 
symmetric periodic orbit with both binary and simultaneous binary collisions in 
the collinear four-body problem. This orbit was proven to exist analytically by 
Ouyang and Yan in [11]. Later, a wider variety of periodic collision-based orbits 
was shown to exist by Shibayama in [12]. The key feature in each of these orbits 
is that they can be reduced to two degrees of freedom. 

In this paper, we study a periodic orbit in the planar rhomboidal four-body 
problem, featuring binary collisions between alternating pairs of masses (as shown 
in Figure 1). This orbit was studied by Waldvogel in [13], who numerically stud- 
ied the chaotic behaviors of the bodies, including escape, quasiperiodic orbits, 
resonances, and quadruple collision. The stability technique of Roberts was used 
in [7] to show that this orbit was linearly stable in the case where all masses are 
equal. Our work extends the stability calculations to the case where there are two 
pairs of equal masses. Under an appropriate rescaling, we need only consider a 
bounded interval of mass values, which simplifies the numerical calculations. 

The remainder of the paper is as follows. In Section 2, we will describe the 
periodic orbit under consideration, and give the coordinate transformation which 
regularizes the binary collisions. Section 3 will outline the standard technique for 
showing linear stability for the orbit, as well as give the result of the symmetry- 
reduction technique developed by Roberts. Section 4 is divided into three portions 
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Figure 1 : The rhomboidal four-body orbit. 



detailing the numerical calculations performed. In particular, the first subsection 
describes how initial conditions for the orbit were obtained, the second gives the 
results for the numerical linear stability computation, and the third describes a 
Poincare section analysis of the orbit. 

2. The Rhomboidal Symmetric-Mass Orbit 

We consider the planar Newtonian 4-body problem with bodies located at 

q\ = (xi,0),q 2 = (0,x 2 ),q3 = ~qiA* = ~12 (2) 

and masses 1, m, 1, m respectively for some m G (0, 1]. The bodies travel along 
the x and y axes, forming the vertices of a rhombus at all times away from colli- 
sion. Binary collisions occur between the bodies with equal masses at the origin. 
For the periodic orbit, the non-colliding bodies have zero momentum at collision 
time. Following collision, the colliding bodies eject along the appropriate coordi- 
nate axis, and the remaining two bodies travel toward collision at the origin, where 
a similar (zero momentum of non-colliding bodies) behavior occurs. (See Figure 
1.) 

Analytic existence of this orbit was shown in [7] in the m = 1 case, and for 



m E (0, 1] in [12]. The Hamiltonian for the system is given by 



1 2 1 2 ' m2 4m 



H = -W\ + - — w 2 — 
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where wi = 2xi,W2 = 2mx2- We can continue the orbits past collision via a reg- 
ularization under which binary collision corresponds to an elastic bounce. To 
regularize these collisions, we use a Levi-Civita-type change of coordinates. Us- 
ing the canonical transformations Qf = xu Pi = 2Q{W{ for i = 1,2, with a change 
of time satisfying dt/ds = x\X2, the regularized Hamiltonian in extended phase 
space is given by 

r = r(Q u Q 2 ,p h P2,E) 



Q\ + Q 2 



This yields the equations of motion 



el = \q\pu 0) 
^ = -^ +e ' m+ (e^-(ei^ +2eiei£ ' (5) 

where ' denotes the derivative with respect to the new time variable s. 

At the time of collision of the two bodies on the x-axis, we have Q\ = and 
P2 = 0. At this time, setting T = yields 
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Hence, Pi = ±8 3//2 , with the sign being the same as the sign on Q2. Similarly, at 
the time of collision of the two bodies on the y-axis, Q2 = P\ = 0. The condition 



r = then gives 



and so P2 = ±(8m 
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The symmetries present in the rhomboidal four-body orbit are the same as 
those presented in the collinear four-body orbit presented in Section 3 of [6]. Let 
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If y(s) is a r -periodic solution to (3) - (6), both Sy(-s + T) and -Sy(-s + T/2) 
are solutions as well. Existence and uniqueness of solutions then imply that 

y(-s + T) =Sy(s) andy(-s + r/2) =-Sy(s) 

for all s. Hence the symmetry group for the rhomboidal four-body orbit is isomor- 
phic to the Klein four group. 



3. Linear Stability 



Note that T is a smooth function defined on IR 4 \ {Q\ =Qi = 0}. Suppose y(s) 
is a r -periodic solution of the system z! = JDT(z), where ' = d/ds, 
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and / and O are the 2x2 identity and zero matrices, respectively. If X(s) is the 
fundamental matrix solution of the linearized equations 



^=jD 2 (y( S m, m 



(8) 



then the monodromy matrix is given by X(T) and satisfies X(s + T) = X(s)X(T) 
for all s. Eigenvalues of the monodromy matrix are also the characteristic mul- 
tipliers of y, and therefore determine the linear stability of y. In particular, y is 
spectrally stable if all of its characteristic multipliers lie on the unit circle, and y 



is linearly stable if it is spectrally stable and semisimple apart from trivial eigen- 
values. 



It is also possible to find the characteristic multipliers of yby using another ini- 
tial condition ^(0) in (8). If Yq is an invertible matrix, let Y(s) be the fundamental 
matrix solution to 

% = JD\y{s)% m=Yo- (9) 

Then by definition of X(s), we have Y(s) = X(s)Y , and so X(T) = Y(T)Y Q l . 
Then we have 

X(T) = Y(T)Y Q l = Y (Y l Y(T)))Y Q l . 

Hence, Yq 1 Y(T) is similar to X(T), so the eigenvalues of Y l Y(T) are the same 
as the characteristic multipliers of y. 



In [3], a technique for factoring the monodromy matrix in terms of symmetries 
of a periodic orbit is presented. These allow the linear stability of a periodic orbit 
to be determined by computing the eigenvalues of the fundamental matrix solution 
at time s = T/N for some natural number N determined by the symmetry group 
of the orbit. Moreover, Roberts in [3] shows that if W is a symplectic matrix such 
that 
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for some matrix K, then the eigenvalues of W lie on the unit circle if and only 
if all the eigenvalues of K are real and lie in the interval [—1,1]. Further, the 
eigenvalues W lie on the unit circle if and only if the eigenvalues of W k lie on the 
unit circle for all k E 7h. Each of these results combine to give the following linear 
stability result. 

Theorem 1. For the rhomboidal symmetric-mass orbit, there exists a choice o/Yq 
such that Yq 1 Y(T) = W 2 for a symplectic matrix W which satisfies (10). More- 
over, K has the form 



K 



-1 
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where Ci represents the ith column ofY{T/A). Therefore, the rhomboidal four- 
body orbit is linearly stable if and only if ^(SJc^) is real and lies within the 
interval [—1, 1], with S as defined in (7). 



The proof of this theorem is precisely sections 3.1 and 3.2 of [6], with greater 
detail presented in [3]. Although the arguments in that section are applied to a 
different orbit, both orbits have the same time -reversing symmetries, and the sym- 
metries of the orbit are the only feature required for the arguments to hold. 



4. Numerical Calculations 

4.1. Determining Initial Conditions 

In order to determine the initial conditions for the rhomboidal orbit, we model 
each of Qi,Qi,Pl, and P 2 by truncated trigonometric polynomials: 

n 

fii = £fl/sin((2/+l)j), (11) 

& = J>sin((2/+l)(j + 7C/2)), (12) 
(=0 

A = i>sin((2/+l)(s-7l/2)), (13) 

(=0 
n 

P 2 = l>sm((2i + 1». (14) 

i=0 

The choice of trigonometric polynomials is natural for modeling periodic behav- 
ior. A similar technique was carried out by Simo in [14]. The time shifts and 
choice of odd-only multiples of s correspond to symmetries of the orbit. In par- 
ticular, for these polynomials, the time-reversing symmetries shown earlier are 
built-in, and the non-colliding bodies have zero net momentum at collision time. 
For a fixed n, we numerically minimize the value of 

J q {(Q[-Q[) 2 + (Q' 2 -Q' 2 ) 2 + (P[-P[) 2 + (P2-Pi) 2 ) ds 

where the minimization is taken over the space of coefficients {a ( -,&;,c;,<i;}. We 
combine this with a root-finding technique to find the appropriate value of E for a 
2n -periodic orbit, as in [5]. Once these trigonometric polynomials are determined, 
we can extract the initial conditions for the orbit by evaluating Q\ , Q 2 , Pi ,P 2 at any 
fixed time s E [0, 2k] . 
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Figure 2: The value of Q2 as a function of m. Figure 3: The value of £ as a function of m. 



We obtain the initial conditions for the 2n -periodic orbit for m = 1 by rescaling 
the conditions given in [7] . It is easy to check that if (Q\ (s) , Q% [s) , Pi is) , P2 is) , E) 
is a solution for (3) - (6), then (eQi(es),eQ2(£-s),Pi(es),P2(£s),E/e 2 ) also satis- 
fies (3) - (6). Moreover, this rescaling does not change the linear stability of the 
orbit. Given the initial conditions in [7] and rescaling with £ w 1.55, the orbit 
is roughly 2n -periodic, which we verify by integration using the standard fourth- 
order Runge-Kutta-Fehlberg algorithm. A standard curve-fitting technique can 
then be used to give the coefficients {a;,£ ; -,c ; -,d,} in the equations (11) - (14). Af- 
ter that, a gradual "step-down" technique (as in [5]) can be used to find the initial 
conditions for other values of m £ (0, 1]. If we assume the initial conditions occur 
at the time the two bodies with mass 1 collide, we have Q\ = P2 = 0. Also, as 
before, we know that P\ = —2 a/2, so we need only to find the values of Q2 and E. 
The results of the numerical calculation are shown in Figures 2-5. 

4.2. Linear Stability Computations 

Using the initial conditions for y(s) leading to a 27T. -periodic orbit given earlier, 
we numerically compute Y(T/4) with the initial condition 
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Figure 5: Orbits in terms of the original physical variables x\ and X2 obtained by coordinate 
transformation of the data shown in Figure 4. 



from Theorem 1 using a standard Runge-Kutta-Fehlberg algorithm with a fixed 
time step 

7l/2 
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The value of c\ (S/C4) plotted against m is given in Figure 6. Numerically, the 
computed values of c\ (S/C4) lie below 1 for all m > .04, and above 1 for m < .03. 
Accounting for numerical error, we have shown the following 

Theorem 2. There exists some positive number £ such that the symmetric mass 
periodic orbit of the regularized planar rhomboidal four-body problem (3) - (6) is 
linearly stable for m G (.03 + £, 1]. 

Note that in the case m = 1, our computation confirms the calculation given in 
[7], wherein an additional symmetry is present in the orbit. 
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Figure 6: The value of c\ (SJC4) as a function of m. 



4.3. Poincare Section Analysis 

To numerically analyze nonlinear stability, we find a suitable Poincare section 
for the orbit. This was done in the m = 1 case in [13]. Our more general Poincare 
section is based on techniques presented in [15] and [9]. For any value of m, we 
seek a number a such that 

— = a 

X2 



is maintained throughout the orbit, with x\ and X2 as defined in (2) earlier. In 
other words, the value of a corresponds to the ratio of x\ and X2 in a homographic 
orbit where the trajectories of the bodies correspond to total collapse (or ejection 
from total collapse). We find this value of a by solving the standard equations 
of motion (1) with the substitutions x\ = OU2 and x\ = ax'2. Doing so, we find 
that the required value of a for a given mass m is the single positive root of the 
12th-degree polynomial 

(1 +oc 2 ) 3 (l -ma 3 ) 2 -64a 6 (m- l) 2 = 0. 

The value of a as a function of m is plotted in Figure 7. 
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Figure 7: The value of a as a function of m. 



For fixed E = — 1 , we define a Poincare section E to be the two-dimensional 
surface given by x\ = OU2 in the phase space defined by the variables (x\,X2,x\,X2). 
(Note that x\ and ±2 are simply linear re-scalings of w\ and W2.) Restricting to 
E = — 1, we find a bound on the possible values of x\ . Specifically, if x\ = X2 = 
on E, the condition E = — 1 requires that 

1 m 2 a Ama 

2 2 y/a 2 + 1 

For a set of initial conditions on E, the requirement E = — 1 necessarily implies 
that x\ < r max , and if either of x\ or ±2 are non-zero, then the strict inequality 
x\ < r max holds. 



We define coordinates (r, 0) on E by 



r = , = tan . 

r ma x \ax 2 J 

Under this change of coordinates, the homographic orbit corresponds to the line 
= tc/4. For a 9 x 15 grid of equally spaced initial conditions in (r, 0) we numeri- 
cally integrate (3) - (6) for the corresponding initial conditions and record the first 
200 intersections of the orbit with E. (Integration was preemptively terminated if 
any of Qi.Pi exceeded 1000 in absolute value.) The results of this are shown in 
Figures 8-12. The observed concentric rings numerically match the predicted 
result of Moser's Invariant Curve Theorem in [16], and show that the rhomboidal 
symmetric-mass orbit is nonlinearly stable for m G (.03 + £, 1] for the same £ as in 
Theorem 2. 
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Figure 8: Poincare sections plotted for m = .1 (top) and m = .2 (bottom). In these plots, r lies on 
the vertical axis. The homographic orbit at 8 = Ji/4 is not plotted for clarity. 
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Figure 9: Poincare sections plotted for m = 



.3 (top) and m = .4 (bottom). 





Figure 10: Poincare sections plotted for m = 



.5 (top) and m = .6 (bottom). 





Figure 1 1 : Poincare sections plotted for m = 



.7 (top) and m = .8 (bottom). 





Figure 12: Poincare sections plotted for m = 



.9 (top) and m = 1 (bottom). 
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